Exploring the geometry of the bifurcation sets in parameter space

By studying a nonlinear model by inspecting a p-dimensional parameter space through \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$(p-1)$$\end{document}(p-1)-dimensional cuts, one can detect changes that are only determined by the geometry of the manifolds that make up the bifurcation set. We refer to these changes as geometric bifurcations. They can be understood within the framework of the theory of singularities for differentiable mappings and, in particular, of the Morse Theory. Working with a three-dimensional parameter space, geometric bifurcations are illustrated in two models of neuron activity: the Hindmarsh–Rose and the FitzHugh–Nagumo systems. Both are fast-slow systems with a small parameter that controls the time scale of a slow variable. Geometric bifurcations are observed on slices corresponding to fixed values of this distinguished small parameter, but they should be of interest to anyone studying bifurcation diagrams in the context of nonlinear phenomena.


Exploring the geometry of the bifurcation sets in parameter space
Roberto Barrio 1* , Santiago Ibáñez 2 & Lucía Pérez 2 By studying a nonlinear model by inspecting a p-dimensional parameter space through (p − 1) -dimensional cuts, one can detect changes that are only determined by the geometry of the manifolds that make up the bifurcation set.We refer to these changes as geometric bifurcations.They can be understood within the framework of the theory of singularities for differentiable mappings and, in particular, of the Morse Theory.Working with a three-dimensional parameter space, geometric bifurcations are illustrated in two models of neuron activity: the Hindmarsh-Rose and the FitzHugh-Nagumo systems.Both are fast-slow systems with a small parameter that controls the time scale of a slow variable.Geometric bifurcations are observed on slices corresponding to fixed values of this distinguished small parameter, but they should be of interest to anyone studying bifurcation diagrams in the context of nonlinear phenomena.
Given a family of dynamical systems, the term bifurcation is used to refer to any qualitative change in the dynamics under variation of parameters.The stratification B of the parameter space by bifurcation points is said the bifurcation set.On the other hand, when one explores a q-dimensional parameter space with (q − s)-parametric families of s-dimensional manifolds S ε , where s < q and ε ∈ R q−s , the geometry of B itself can lead to changes in the sets B ε = B ∩ S ε .We are interested in these changes, that we call geometric bifurcations.As might be expected, the study of geometric bifurcations is closely related to the Theory of Singularities.
In this article we use as test examples two classical mathematical neuron models: the Hindmarsh-Rose and the FitzHugh-Nagumo systems.The main reason is that these models, as examples of fast-slow systems, have distinguished parameters that dictate a specific way of exploring the parameter space.Furthermore, in multiscale systems there are bifurcation sets in the parameter space where different structures are exponentially close to each other 1,2 , so there are more types of geometric bifurcations, as will be described later.But we emphasize that they can be present in any parametric exploration of a system regardless of whether it is fast-slow or not.
The Hindmarsh-Rose 3 (HR) system was introduced as a reduction of the Hodgkin-Huxley equations 4 for the propagation of nerve impulses in axons.The HR model is described by three coupled nonlinear ODEs: where x is the membrane potential, y the fast and z the slow gating variables for ionic current.Typically, b, I and ε are considered as free parameters, whereas the remaining ones are set as follow: a = 1 , c = 1 , d = 5 , s = 4 , x 0 = −1.6.
The FitzHugh-Nagumo system 5,6 is another simplified version of the Hodgkin-Huxley model 4 , and the reduced ODE system can be written in the form 7 : (1) To illustrate the main aim of this article, the spike-counting (SC) technique 8,9 is used in Fig. 1 to show regions of periodic tonic spiking, chaotic bursting, and regular bursting with different number of spikes for the classical Hindmarsh-Rose neuron model (more details will be provided later).The SC sweeping method calculates for each set of parameter values and initial conditions the number of loops (spikes) of the stable limit cycle (when it exists), and this number is color-coded according to the number of spikes (vertical bar).Specifically, the dark blue regions correspond to stable spikings.The strips of a spike-adding cascade, ranging from blue to green, are related to 10 fold/Hopf and fold/hom bursting, which becomes chaotic 11 in a chain of onion-like bulbs depicted in brown.Additionally, several bifurcation lines are superimposed, such as the black one for homoclinic bifurcations and the red curves for period-doublings of periodic orbits.It is worth noting that the bifurcation curves of periodic orbits arise from homoclinic bifurcation points of codimension two, namely, the HR model exhibits 1 orbit-flips (OF), inclination-flips (IF) and Belyakov bifurcations.Furthermore, we highlight how the structure of the parameter plane (b, I) changes as the small parameter ε increases, moving away the system from the singular limit at ε = 0 .An example would be the disappearance of some color bands and several bifurcation points of codimension two (e.g., the green and pink points) while, on the contrary, connections appear between bifurcation curves of periodic orbits (in red).Finally, we also look at how the Z-shape of the homoclinic curve (in black) evolves and disappears.Regarding these phenomena, the question that naturally arises is what is the mechanism underlying all these changes.As we will see later, the answer is that the bifurcation set goes through several geometric bifurcations as ε varies.
The relevance of geometric bifurcations has already been pointed out in the literature, although the terminology we are using in this paper is new.Geometric bifurcations explain the creation of isolas of homoclinic and heteroclinic bifurcations.Namely, Algaba et al. 12 study Chua's equation and show how exploring the 3-parameter space with planar slices, a pair of codimension-two bifurcation points (T-points) collapse and disappear.When the collision occurs, the intersection between the slice and the codimension-two bifurcation curve is non-transversal.Through this process, spirals of codimension-one bifurcations approach together, giving rise to isolas.The same authors provide 13 a theoretical model explaining the behavior.They also study 14 the distribution of homoclinic bifurcation curves in the Lorenz system, and illustrate again how the evolution of a bifurcation set depends on the way in which the parameter space is explored.Wieczorek and Krauskopf 15,16 study a threeparameter family of three-dimensional vector fields that models an optically injected laser.Planar slices provide images of the bifurcation set where codimension-two bifurcation points collapse together, and codimension-one bifurcation curves contact each other.Among other phenomena, isolas of homoclinic bifurcations are exhibited.The systems considered in the papers we have just mentioned [12][13][14][15][16] and the neuron models we will consider later give an idea of the types of scenario in which geometric bifurcations can occur.In general, they can be useful in any model that depends on three (or more) parameters when the parameter space is explored using 2-dimensional slices.Therefore, the notion of geometric bifurcation should be of interest to anyone studying bifurcation diagrams in the context of nonlinear phenomena.
In this paper, we provide a minimal theoretical framework regarding geometric bifurcations by appealing to the theory of singularities on differentiable manifolds, and particularly, to Morse Theory 17,18 .More importantly, it is shown how all the bifurcations that we introduce are observed in rather common neuron models.As already mentioned, we provide illustrations for the Hindmarsh-Rose 3 and FitzHugh-Nagumo systems 5,6 .Both are fastslow systems obtained as simplifications of the Hodgkin-Huxley model 4 for the propagation of nerve impulses in axons.In both cases there are very rich bifurcation sets that, for instance, include homoclinic bifurcations of codimension one and two.In addition, we will see how geometric bifurcations arise when the parameter space is explored with slices corresponding to fixed values of a distinguished parameter ε , namely, the parameter that controls the time scale of the slow variable.These models illustrate a key fact, the study of geometric bifurcations only makes sense when working with families where there are distinguished parameters that determine a specific way of inspecting the parameter space; otherwise, either a reparametrization may remove singular points or tangencies can be avoided by choosing a different direction for sweepings.
Golubitsky and Schaeffer 19 also use the term "distinguished parameter" to refer to a specific parameter whose variation can be seen as typical or natural in a given model.They work with scalar equations of the form for a single unknown x and with α ∈ R and β = (β 1 , . . ., β q−1 ) ∈ R q−1 .Typically, (3) stands for the steady-state equation of a certain dynamical problem.In this setting, they refer to α as a distinguished parameter and say that β 1 , . . ., β q−1 are unfolding parameters.Golubitsky and Schaeffer wonder about how the sets change as β varies.In this context, they characterize and classify the bifurcation points (x 0 , α 0 , β 0 ) -singularities in their terminology-where a change occurs.The notion of geometric bifurcation is also close to this theory of singularities.We consider families of dynamical systems depending on q parameters ( , ε) ∈ R q−s × R s , with s < q and explore the parameter space taking (q − s)-dimensional manifolds with ε fixed.We wonder about how the intersection between these manifolds and the bifurcation set changes as ε varies.In our setting, there is no state variable and ε represents the distinguished parameter (or parameters in case that s > 1 ).For the slow- fast models that we consider in this paper, the choice of distinguished parameter is indeed a "natural" one: the small parameter ε that controls the behavior of the slow variable.Note that the theory of singularities 19 allows to understand, for instance, the appearance of isola of equilibria in bifurcation diagrams 20 .
The notion of geometric bifurcation incorporates a rather particular concept of codimension, which we borrow from Wieczorek and Krauskopf 15,16 .Thus, bearing in mind a three-dimensional parameter space, we can find only three types of geometric bifurcations regarding codimension: 1. Codimension one-plus-one geometric bifurcations: Those that occur on a surface of codimension-one bifurcations when the parameter space is explored with one-parametric families of two-dimensional manifolds.2. Codimension two-plus-one geometric bifurcations: Those that occur on a curve of codimension-two bifurcations when the parameter space is explored with one-parametric families of two-dimensional manifolds.3. Codimension one-plus-two geometric bifurcations: Those that occur on a surface of codimension-one bifurcations when the parameter space is explored with two-parametric families of one-dimensional manifolds.
The structure of the paper is as follows.First, simple pictures that allow one to illustrate the notion of geometric bifurcation are provided.They help to understand the particular way in which the codimension of a geometric bifurcation is specified.The existence of a codimension two-plus-one geometric bifurcation has subsidiary consequences on the codimension-one bifurcation sets that emerge from it, and these effects are also discussed.Later, we formalize some notions and provide a minimal theoretical support for dealing with geometric bifurcations.Specifically, codimension one-plus-one and two-plus-one bifurcations are discussed.After that, we show how geometric bifurcations arise in the aforementioned neural models.The provided preliminary illustrations and subsequent theoretical schemes should be useful to understand many of the phenomena exhibited in the neuron models.Finally, some conclusions are provided.

Geometric bifurcations: illustrative examples
Let X 1 , 2 ,ε be a three-parameter family of vector fields, and assume that a bifurcation surface When we explore the parameter space taking horizontal planes with ε fixed, we observe that, for ε < 0 , there is a bifurcation curve, the circle 2 1 + 2 2 = −ε (see top panel in Fig. 2).In the absence of additional bifurcations, there is one equivalence class for values of the parameters inside the circle and one outside.When ε = 0 , the interior class disappears because the circle collapses to a point.Thus, when ε > 0 , bifurcations are no longer observed.If the bifurcation surface is given by , we observe a different behavior (see bottom panel in Fig. 2).For ε < 0 , there are two bifurcation curves, (the two branches of the hyperboloid ε = 2 1 − 2 2 ).In the absence of additional bifurcations, there exist two equivalence classes.When ε = 0 , the bifurcation set is given by 2 1 − 2 2 = 0 , two straight lines, and there is a singular point at (0, 0).For ε > 0 , we find two curves again, but placed on differ- ent sectors.Note that in the passage through ε = 0 , the class which corresponds to the neighborhood of (0, 0) changes.In both cases, we say that there is a one-plus-one geometric bifurcation at ε = 0 .Note that intersections between the horizontal planes and the surface M are all transverse except for ε = 0.
In the previous examples, a geometric bifurcation appears because the three-parameter space is explored with planes and there is one which is tangent to a bifurcation surface.If there exists a codimension-two bifurcation curve exhibiting a point of tangency with one of planes, we say that there is a codimension-two-plus-one geometric bifurcation.Both panels in Fig. 3 show a codimension-two bifurcation curve (blue) C. When the parameter space is explored with horizontal planes (fixed values of ε ), there is a value of ε for which the plane has a quadratic tangency with C at a point p.Two codimension-two bifurcation points are exhibited for lower values of ε and no one for bigger values.We say that there is a codimension-two-plus-one bifurcation at p.
More important, another interesting question is to wonder about the changes that occur in the bifurcation strata that arise from other strata of lower codimension.Both panels in Fig. 3 show a surface M (grey color) of codimension-one bifurcation points with a boundary given by the curve C. In the left panel, where the point p is a saddle point on the surface, we can see how, as ε increases, the horizontal slides show two curves (red color) that join together at p to form a unique curve.On the other hand, looking at right panel we see that the point p is  www.nature.com/scientificreports/ a maximum on the bifurcation surface.In this case, as ε increases, the horizontal slides show a unique curve (red color) that collapses at p and disappears.The changes on the red curve are a subsidiary effect of the codimension two-plus-one geometric bifurcation.Note that in both cases, the bifurcation curves that are shown in each slide emerge from the codimension-two point with a well-defined tangent and we can observe how, "generically", they reconnect with each other or simply disappear, when the codimension-two points are no longer present.The former case is the one exhibited by the red curves in Fig. 1.The coupled laser model studied by Blackbeard 21 also provides similar reconnections (see Figures 4.6 and 4.7 in his thesis).The case of spirals of bifurcations emerging at codimension-two bifurcation points has also been studied 12,13,15,16 New scenarios arise if we reduce the dimension of the slices.Let X 1 ,ε 1 ,ε 2 be a three-parameter family of vector fields and assume that Surface M is shown in Fig. 4a together with a horizontal projection plane (green colors), that is, a plane corresponding to parameters (ε 1 , ε 2 ) .On this two-dimensional parameter plane, we describe the unfolding of the geometric bifurcation.Each time we fix a pair (ε 1 , ε 2 ) , we can consider a vertical line that will meet M at either, one, two or three points.All the possible behaviors are shown in Fig. 4b, although lines are depicted horizontal.If (ε 1 , ε 2 ) belongs to the interior of region 3, i.e., the dark-green region bounded by the curve 27ε 2 1 + 4ε 3 2 = 0 (in red), we find three bifurcation points a, a ′ and a ′′ (Fig. 4b-line 3).When parameters (ε 1 , ε 2 ) cross line 2 (Fig. 4b), a and a ′ collapse at a codi- mension-one-plus-one bifurcation displayed in the two-parameter space we obtain by fixing ε 2 (see point b).The same occurs along line 2 ′ , but the collision is between a ′ and a ′′ instead (see point c ′ ).When parameters (ε 1 , ε 2 ) are in region 1 (Fig. 4b-line 1), we find only one bifurcation point at d, points b and c ′ disappear when param- eters cross lines 2 and 2 ′ , respectively.We say that there is a codimension-one-plus-two geometric bifurcation at (ε 1 , ε 2 ) = (0, 0) .Note that Fig. 4b also shows the case corresponding to that bifurcation point itself (the pink point).The Z-shape, exhibited by the bifurcation curves given by the intersection between M and the vertical planes with ε 2 < 0 fixed, should be compared with the evolution of the black curves in Fig. 1.
Remark 1 In principle, any dynamical behavior can disappear when traversing a geometric bifurcation.Consider, for instance, Figure 2 (top), where the dynamics within the light blue region may correspond to the presence of a specific attractor, potentially within a framework of multistability.However, this dynamic stops happening for slices where ε > 0.

Geometric bifurcations and Morse theory
Recall that a codimension-k bifurcation is generic in the set of families of vector fields dependent on a number of parameters larger than or equal to k, but it can be avoided in families dependent on fewer parameters.A codimension-k bifurcation is characterized by a set of k degeneracy conditions independent of each other and a set of non-degeneracy (open) conditions.Any family satisfying these conditions is said to be a generic family or a generic unfolding.
Let X ,ε : R l → R l be a C ∞ family of vector fields with ( , ε) ∈ R q × R .Assume that the family is a generic unfolding of a codimension-k bifurcation at ( , ε) = (0, 0) , with k ≤ q .Let M be the (q + 1 − k)-dimensional smooth manifold through (0, 0) where such bifurcation is exhibited.Definition 2 We say that the bifurcation point (0, 0) in the smooth manifold M is geometrically generic with respect to ε if the hyperplane ε = 0 is transversal to M at (0, 0).Otherwise, we say that the bifurcation point is geometrically degenerate with respect to ε.
Let us recall the notion of Morse function 22 .Definition 3 Given a smooth manifold M and a smooth function f : M → R , we say that a critical point p 0 of f is nondegenerate if its Hessian is nondegenerate.We say that f is a Morse function if all critical points are nondegenerate.

Definition 4
We say that the geometrically degenerate bifurcation at (0, 0) has geometric codimension-one with respect to ε if the height function is locally a Morse function.
The condition for geometric degeneracy is the fact that (0, 0) is a critical point of the height function P ε .Imposing that P ε is a Morse function, we are setting a generic geometric condition.Following Wieczorek and Krauskopf 15,16 , we introduce the notion below.

Definition 5
We say that a bifurcation point has codimension n-plus-1 if it corresponds to a bifurcation of codimension-n from a dynamical point of view, but the same point has codimension-1 from a geometrical point of view.

Remark 6
To provide a theoretical framework for geometric bifurcations with codimension m, with m > 1 , one needs to consider m distinguished parameters (ε 1 , . . ., ε m ) .In that context, the notion of codimension n-plus-m bifurcation would make sense.Nevertheless, in this paper, we do not enter into this formalization.However, note that in the previous section we have illustrated the case of a codimension-one-plus-two geometric bifurcation.

Remark 7
For generalisations, one should go a step further with respect to classical Morse Theory and refer to Cerf Theory 23,24 on families of real-valued functions.However, since we want this notion of geometric bifurcation to reach a wide and heterogeneous audience, a deeper investigation at the theoretical level is out of the scope of this paper.
In the sequel we assume that q = 2 .Suppose first that there exists a codimension 1+1 bifurcation point at (0, 0) ∈ R 2 × R .Since the height function is a Morse function locally around (0, 0), it follows from the Morse Lemma that using convenient coordinates = ( 1 , 2 ) , one can write either ε In the first two cases, the level curves are circles (isolas) for ε > 0 (resp.ε < 0 ) and empty sets for ε < 0 (resp.ε > 0 ) (see Fig. 2(top)).We say isola-type to refer to this bifurcation.In the latter case the level curves correspond to a saddle case (see Fig. 2(bottom)) and we say simply saddle-type bifurcation.These two cases are also discussed by Wieczorek and Krauskopf 15,16 .The isola-type and saddle-type geometric bifurcations correspond to the codimension-one singularities named isola-center and simple bifurcation, respectively, by Golubitsky and Schaeffer 19 .Nevertheless, as we already mention in the introduction, there exist singularities which do not match any geometric bifurcation, the hysteresis point 19 , a codimension-one singularity, is one example.Any two isola-type geometric bifurcations are equivalent in the sense that, using appropriate coordinates, their respective height functions can be written in a unique canonical form, either with a maximum or a minimum.The same can be said for saddle-type geometric bifurcations.Of course, we are assuming that the underlying dynamical bifurcation is the same.
Regarding codimension-two-plus-one bifurcations in three-parameter spaces, in principle, we only distinguish one type.Indeed, given a curve C corresponding to a codimension-two bifurcation with a folding point, we can apply again the Morse Lemma and choose a convenient coordinate µ along M such that ε = ±µ 2 .There- fore, the level sets are given by two points that collapse and disappear (compare with the illustration provided in Fig. 3).However, when we consider the bifurcation surfaces of codimension one that may emerge from C , different types can appear (see again Fig. 3).Note that the typology can include two cases or more, since there can be two or more bifurcation surfaces emerging from C , even an infinite number.
Another relevant point is to study how the different geometric bifurcations can appear in the bifurcation diagram of a given family of dynamical systems.We need to understand how a bifurcation surface exhibiting geometric bifurcations can look like.Figure 5 shows a scheme of different scenarios.A key point is the topology of the bifurcation surface, having one or several tubular structures or even maxima.We classify the different 2D manifolds using the Morse theory 18,22 .In all cases the value of a parameter ε on the surface is our Morse function.
The different possibilities for Morse manifolds shown in Fig. 5 are illustrated in next section, namely, the Hindmarsh-Rose (resp.Fitzhugh-Nagumo) model exhibits Cases I to III (resp.Case IIb).-Case I-surface is topologically equivalent to a cylinder, that is, a two-holes sphere.All values of ε are regular and hence the struc- ture of the surface is the same, a single circle, for all level sets.-Case II-surface, a three-holes sphere, is called a "pair of pants" surface in the context of topology.In this case the Morse function given by the value of ε has a saddle point (the red point in the middle plot of Fig. 5).This point is a codimension-one-plus-one geometric www.nature.com/scientificreports/bifurcation point.In this case, the passage through the saddle corresponds to a connecting transition: two isolas collide in the saddle point and give rise to a unique isola.In the case shown in-Case III-, the pair of pants structure is broken, giving rise to two independent surfaces that disappear finally at a maximum of the surface when ε grows (see the blue points).This case is a two discs shape (two one-hole spheres).Again the points of maxima are codimension-one-plus-one bifurcation points.Finally,-Case IIb-, that we call "duck-foot surface", is obtained from the combination of a "pair of pants" and a "disc" (two Morse surfaces that meet along a curve).Note that this case requires the existence of a codimension-two bifurcation curve (magenta color) that acts as a limit of the tubular structures and enables the continuation in a disc surface.Therefore a cut on a transversal plane gives us the picture of two isolas connected with a curve.More combinations of Morse surfaces are possible, but we only discuss the one detected in the Fitzhugh-Nagumo system.Bifurcation surfaces include codimension-two bifurcation curves (black and blue lines) and codimensionthree bifurcation points (black points).The value of ε along the curves of codimension-two bifurcation plays a crucial role in the understanding of the global picture.Therefore, this parameter has a relevant "physical" meaning as it is assumed to be the parameter that can be changed in the system giving different phase space dynamics.All points along a curve of codimension-two bifurcation where ε reaches an extremum (white and green points) are codimension-two-plus-one geometric bifurcations.

CASE I C ASE II CASE III
cylinder-shape pair of pants-shape 2 x discs-shape (2-holes sphere) (3-holes sphere) ( 2 x 1-hole sphere) codimension (2+1) -"visible fold" codimension (2+1) -"invisible fold" codimension-three bifurcation point codimension (1+1) -"saddle-type" codimension (1+1) -" isola-type" surface of codimension-one bifurcations curves of codimension-two bifurcations codimension (2+1) -"saddle-type"+ "visible fold" www.nature.com/scientificreports/ The tubular structures exhibited by the Morse surfaces allow us to distinguish between two different types of bifurcation curves of codimension two: those that appear in pairs, on opposite sides of the cylinder; and those that appear as single curves with a fold point where they cross from one side of the cylinder to the other.All curves shown in Case I and those in blue color that we observed in the pair of pants of cases II and IIb are of the first type.The second type corresponds to the bifurcation curves in black color shown in Case II and all the curves in Case III.An independent case is the bifurcation curve (blue color) which is contained in the disk component of the Morse surface in Case IIb.In all cases, the points on the curves where ε reaches a maximum on the curves are codimension-two-plus-one geometric bifurcations.White points are folds with both branches contained in a unique side of the cylinder (or inside the disk for the Case IIb).Green points are folds where the two branches belong to different sides of the cylinder.

Remark 8
It is important to point out that along the red curve the-duck foot is not a regular surface.In particular, Definition 4 does not apply to the red-white point shown in the figure.However, as already mentioned, the duck-foot surface must be understood as two Morse surfaces glued together along a curve.Looking at the pair of pants, the point corresponds to a codimension-one-plus-one saddle-type geometric bifurcation (which explains the red color).But if one thinks in the disk, or better in the red curve, the point corresponds to a codimension-two-plus-one geometric bifurcation (which explains the white color).
Let d be the distance between both sides of a tubular structure, as depicted in Fig. 5.If d is small enough, the two branches arising from green folds are indistinguishable and the folding point may be wrongly perceived as an end point.This situation happens in the fast-slow examples provided in the next section, where d is exponentially small.To distinguish these "false" end points, green points are said to be "invisible folds", in contrast with white points that are said to be "visible folds" (note that this situation is for d ≪ 1 ).In the next section, we do not include models with d big, but they would not provide any other noteworthy behavior either.Note that the appearance of invisible folds is linked to the formation of the "pants" or the splitting of a Morse surface in different components.This allows two disconnected isolas of a codimension-two bifurcation curve to meet and create a unique isola.Also, note that there are infinite options to create different geometric surfaces (not regular surfaces) mixing different Morse surfaces and codimension-two bifurcation curves.For example, in Fig. 6 we show some additional possibilities.In case (a) we illustrate the case where the bifurcation surface is broken due to different codimension-two curves (this situation can occur, for example, in Hopf bifurcation surfaces broken by Bogdanov-Takens bifurcation lines).Case (b) presents the "pair of pants" ended in two maxima.Case (c) shows two "discs" connected by a cylinder glued together by codimension-two curves.Case (d) is built by adding more and more "pair of pants" with different configurations.Note that all the structures can be repeated as many times as we want.The final example (e) shows the case of having a codimension-2 curve that gives rise to a countable pencil of secondary bifurcations that split from the bifurcation line "as pages of a book" (for example, in orbit-flip or inclination-flip bifurcations of type C 25 where countable pencils of period-doubling and fold bifurcations of periodic orbits occur).In all these combinations geometric bifurcations of those described above may appear.

Some models exhibiting geometric bifurcations
In this section we will consider two models that exhibit geometric bifurcations, namely, the Hindmarsh-Rose model and the FitzHugh-Nagumo system, showing how geometric bifurcations permit visualizing the global parameter space.
As discussed in the introduction, the main reason for using these models is that fast-slow systems naturally have distinguished parameters that dictate a specific way of exploring the parameter space, and there are sets of bifurcations in the parameter space where different structures are exponentially close, but note that geometric bifurcations may be present in a parametric exploration of any dynamical system.

The Hindmarsh-Rose model
The Hindmarsh-Rose model (1) has been deeply studied in recent years 1,2,8-11,26-28 using different techniques.The system is particularly useful to understand the bursting phenomena and, in particular, the mechanisms of spike-adding.Behind these dynamics there are homoclinic bifurcations 8,9,11,27,28 .The global homoclinic structure and also part of the tangle of bifurcations of periodic orbits that emerge from the homoclinic skeleton have been unveiled 1,26 , and this allows to explain the spike-adding processes 1, 10 .
In most of fast-slow systems with explicit small parameters, these parameters play a significant role and under their variation drastic changes in the global phase space are exhibited.The HR model is a paradigmatic example of this fact.As the small parameter ε increases, numerous changes occur.In what follows, we show how the HR model exhibits geometric bifurcations with respect to ε that, linked with dynamic bifurcations, explain these changes.Therefore, the HR model shows how geometric bifurcations help us visualize the global bifurcation diagram in the parameter space.
From previous analysis 1 , it follows that there exist many codimension-one homoclinic bifurcation surfaces which are exponentially close to each other and its number grows to infinity when the small parameter tends to zero.Moreover, since these homoclinic surfaces are tubular, the intersection of each surface with horizontal planes produce isolas (closed curves).These isolas exhibit a pair of extremely sharp folds and their width is also exponentially small (compare with the scheme given in Fig. 5 with d ≪ 1 ).Folding points determine two differ- ent sides in the isola and also in the bifurcation surface.Typically, for parameter values on one of the "sides" of the homoclinic bifurcation isola, the homoclinic orbit exhibits n spikes and, for parameter values on the another face, n + 1 .This explains the notation hom (n,n+1) to refer the different isolas and surfaces.In Fig. 7 we show three-parameter plots showing some codimension-one homoclinic bifurcation surfaces and codimension-two homoclinic bifurcation curves computed using the continuation AUTO software 29,30 .To construct the surfaces, a collection of codimension-one homoclinic bifurcation curves are computed for different values of ε .In plots (a), (b) and (c) the hom (1,2) , hom (2,3) and hom (11,12) surfaces are given, respectively, which include Belyakov and Inclination-Flip codimension-two bifurcation curves.This partial bifurcation diagram illustrate, respectively, Cases I, II and III shown in the theoretical scheme provided in Fig. 5.
Looking for codimension-one-plus-one geometric bifurcations, we pay attention to the bifurcation surfaces themselves.In Fig. 7c we see how hom (11,12) splits into two disconnected components and each of them has a maximum (with respect to ε ).These two points are isola-type codimension-one-plus-one geometric bifurcations with respect to ε .Also, as shown in plot (b) for the surface hom (2,3) , a saddle-type codimension-one-plus-one geometric bifurcation is detected on the bifurcation surface.Most likely, isola-type codimension-one-plus-one geometric bifurcations are present in all homoclinic surfaces if the small parameter ε grows enough and it is no longer "small".This fact explains why, as the small parameter grows, fewer bands of color appear in the 2D plots of Fig. 1 for the largest value of ε , indicating that the bursting orbits have fewer spikes.
On the other hand, codimension-two-plus-one geometric bifurcations correspond to folds, with respect to ε , that appear along curves of codimension-two homoclinic bifurcations and, as explained in the previous section, we distinguish between "visible" and "invisible" folds.Recall that in all two-dimensional manifolds of codimension-one homoclinic bifurcations we distinguish two leaves which are exponentially close (in fact the two leaves that form the isolas) and therefore they are indistinguishable in our visualization of the numerical results.hom (2,3) hom (11,12) invisible folding curve isola of Belyakov bifurcations  2) , hom (2,3) and hom (11,12) , respectively.Some codimension-two homoclinic bifurcation curves are shown over the surfaces and some geometric bifurcations are pointed.
They glue together along curves of sharp folding marked with a red line in Fig. 7.When a curve of codimensiontwo homoclinic bifurcation folds inside one of the leaves, we get a codimension-two-plus-one visible fold, but, if the fold is along one of the curves of folding of the whole surface (that is, going from one leaf to the another), and hence it is hidden for visualization, we have a codimension-two-plus-one invisible fold.Note that visible folds appear in pairs, one on each leaf, but invisible folds are unique points.On the surface hom (1,2) , both the Belyakov and the Inclination-Flip bifurcation curves show "visible" folds.On the surface hom (2,3) , the Belyakov bifurcation curve also undergoes a "visible" fold, but the Inclination-Flip curve presents now an "invisible" fold.And on the surface hom (11,12) , both curves exhibit "invisible" folds.The phenomenon of the "invisible" folds is a direct consequence of the existence of very thin tubular structures (and therefore isolas when a two-parameter section is considered) where the two leaves are infinitesimally close each other.Thus, if a codimension-two homoclinic bifurcation curve reaches the homoclinic surface folding curve, then it continues to the other side because the conditions in the phase space that are required to have the codimension-two bifurcation are still satisfied on the other side.
In the Introduction, we show in Fig. 1 how some codimension-two points disappear, and later we also explain that a subsidiary effect is the fact that different curves of periodic bifurcations may connect together, as illustrated in Fig. 1c.Now we can connect all these phenomena with some geometric bifurcations.In Fig. 8 we show some codimension-two-plus-one bifurcations in the HR model increasing the parameter ε .We illustrate schematically the global structure before and after the codimension-two-plus-one bifurcations and how the size of the homoclinic isolas determines the type: "visible" or "invisible" fold.As ε increases, the first geometric bifurcation in Fig. 1 is a codimension-two-plus-one bifurcation that occurs on the Inclination-Flip (IF) homoclinic bifurcation curve.For the hom (1,2) case, as the homoclinic surface is big enough, there is a maximum of each of the IF curves (one on each leaf of the tubular surface), and so we have a "visible" fold and we observe the geometric "collision" of two pairs of IF points, one pair on each of the leaves of the homoclinic surface (the white points on the Figs. 5 and 7).On the contrary, for the rest of homoclinic surfaces, the IF bifurcation curves have "invisible" folds, as they are smaller and the surfaces are composed of isolas disconnected for small values of the parameter.Moreover, when n grows enough hom (n,n+1) has only one branch of the IF curve, that is, IF points are only present on one of the two disconnected components of the isolas for a small value of ε .In these cases, the .Schematic processes of some codimension-two-plus-one bifurcations in the HR model increasing the parameter ε .On the top, SC pictures with some bifurcations illustrate the global panorama at several values of ε .On bottom pictures some numerical bifurcation curves provide an outline of different ways, through "visible" vs. "invisible" folds, to arrive to similar global situations.Note that ε increases from left to right.
IF point seems to disappear on the limit of a homoclinic curve, and what really happens is that it geometrically "collides" with the corresponding point of the another leaf of the surface that we cannot see (the green points on the Figs. 5 and 7).The evolution of the Belyakov curve is similar but for higher values of the small parameter ε. Figure 8 also illustrates how codimension-one bifurcation curves arising from codimension-two bifurcation points are affected by the presence of a geometric bifurcation.For instance, following the surface hom (1,2) , first a collision of two IF points is observed (a codimension-two-plus-one geometric bifurcation) which is accompanied by the reconnection of two branches of period doubling bifurcation curves.
Finally, to complete the panorama of geometric bifurcations, we show in Fig. 9 how the HR model exhibits a codimension-one-plus-two geometric bifurcation point (compare with Fig. 4, where the ε parameter of HR is the parameter ε 2 ).Recall that, in spite of the notion of geometric codimension two has not been formally introduced, we have linked that concept to conceiving the parameter space as a two-parameter family of oneparameter bifurcation diagrams.So, when we describe this codimension-one-plus-two geometric bifurcation point, geometric bifurcations must be understood as tangencies or degenerated tangencies between straight lines in the parameter space (where the value of two parameters is fixed) and a bifurcation surface.In Fig. 9a we show the projection of the hom (1,2) bifurcation curve for different values of ε on the (b, I) parameter plane.The loca- tion of the codimension-one-plus-one bifurcation points ("visible folds") on this projection plane is given and we see how they form a pair of geometric bifurcation curves.These curves have a cusp-type contact giving rise to a codimension-one-plus-two geometric bifurcation point.A scheme of the geometry of the hom (1,2) bifurca- tion surface is given in the plots shown in Fig. 9b presenting a cusp catastrophe geometry giving rise to this codimension-one-plus-two bifurcation point (compare with the illustration provided in Fig. 4).Note that as the curves of homoclinic bifurcation finish cutting this codimension-one-plus-one curve, the cusp catastrophe surface is incomplete in some parametric regions giving just one fold point ("C" shape vs. "Z" shape).Plots in Fig. 9(c) show SC pictures with the hom (1,2) bifurcation curves for three values of ε to see more clearly the process, and illustrating the changes in the geometry of the homoclinic curves.

The FitzHugh-Nagumo system
The FitzHugh-Nagumo system (2) has been extensively studied in the literature 7,[31][32][33][34][35] .Taking s and p as bifurcation parameters, it was shown 7 the existence of C-shaped curves of homoclinic bifurcations (travelling waves in the original PDE FitzHugh-Nagumo system correspond to homoclinic orbits in system (2)).In fact it was observed that these C-shaped curves were isolas, that authors called "homoclinic bananas", of an exponentially small wide d, that is, exhibiting sharp turning points and being d ≪ 1 (see also the explanations regarding invisible folds in Fig. 5).They also showed how "bananas" could split due to the existence of codimension-two bifurcation points on the homoclinic curve that played the role of terminal points for the branches following after the sharp turning points.That is, a "banana split" consists of two isolas joined by a curve.Geometrical explanations for the sharp turns has been provided 31,33,34 and the transition has been analytically described using geometric singular perturbation theory and blow-up techniques 32 .
We will show that the bifurcation diagram of system (2) with respect to parameters (α, s, ε) exhibits geometric bifurcations which are related with the split of the homoclinic banana.
In Fig. 10 we present a three-parameter plot of the codimension-one homoclinic hom (1,2) bifurcation surfaces for the system (2) with p, γ = 0 , = 1 , and considering α , s and ε as bifurcation parameters.Plot (a) shows the theoretical homoclinic structure that corresponds to -Case IIb-, "duck-foot surface", which is obtained from the combination of a "pair of pants" and a "disc" (see Fig. 5).In this case a codimension-two Belyakov bifurcation curve acts as the limit bound of the tubular structures with a secondary homoclinic (after obtaining an extra spike) that returns to the Belyakov points and enables the continuation in a disc surface.That is, a transversal cut of the structure gives an isola for large values of ε ; for an intermediate value of ε , it gives a point of generation of two coupled isolas that coincides with the geometric maxima of the Belyakov curve (a codimension-two-plusone bifurcation); and for small values of ε , it gives a set of two connected isolas, the "homoclinic banana split" that is described in the homoclinic literature (see also the description provided in Remark 8).Plots (b) and (c) provide numerical results obtained using the continuation software AUTO, and they show the homoclinic bifurcation curves for several values of the small parameter ε and the codimension-two Belyakov bifurcation curve.Plot (b) is given in the three-parameter space (α, s, ε) and apparently we only observe just one curve for each ε , and no loop at all, but this is due to the very small distance between the sides of the isolas ( d ≪ 1 ).In order to observe the isolas and the connected isolas, we have to use the AUTO L 2 -norm as it is shown in plot (c) using the parameters (α, ε) and the L 2 -norm.

Discussion
In this paper, partially bringing together previous studies, we introduce the concept of geometric bifurcation and illustrate how it appears in the prominent context of neural models.Of course, the Hindmarsh-Rose and FitzHugh-Nagumo systems are not unique examples.Geometric bifurcations are also expected to appear in other neural models and also in models coming from contexts other than Neuroscience (see examples in the already mentioned references [12][13][14][15][16] ).As we have already mentioned, the contexts in which one can expect the appearance of geometric bifurcations are diverse.In the case of 3-parametric spaces explored with 2-dimensional slides (a scenario ubiquitous in the literature), we could find geometric bifurcations of codimension 1+1 and 2+1.In this paper, geometric bifurcations are introduced in a very general way.The notions are given in a unified manner, for arbitrary families of vector fields, and the most basic concepts and results are formally introduced at the level of geometric bifurcations of codimension one.The two neuron models we consider are of great relevance and for this reason they have been analysed in numerous papers.However, no discussion is available in the literature about geometric bifurcations that appear in these models when a three-parameter space is explored.We provide a description of these occurrences.
It should be noted that the study of geometric bifurcations may be of special interest in the framework of fast-slow systems, but we remark that they can be present in a parametric exploration of a dynamical system regardless of whether it is fast-slow or not.In this case, those parameters that modulate the slow dynamics are distinguished parameters.These parameters are the ones which one should fix get 2D slices.Therefore, they specify a very concrete way of exploring the parameter space.In fact, it is possible to observe variations in the bifurcation diagrams when the distinguished parameters change, that is, signs of the presence of geometric bifurcations.In summary, if a multi-parameter space needs to be explored, there is no other option than to work with k-parameter slices with k = 1, 2 or, at most, k = 3 , and the different phenomena that we explain in this article could surely be present.
We formally introduce the notion of geometric bifurcation.In particular, we utilize a terminology (codimension n-plus-m) that emphasizes the distinction between bifurcations in families of dynamical systems and in "families of families".Moreover, we studied the geometric bifurcations of codimension-one-plus-one and two-plus-one in three-parameter spaces.Additionally, we described a codimension-one-plus-two geometric bifurcation.We were able to explain the appearance of different geometric bifurcations by the combined use of the Morse classification of 2D manifolds and bifurcation theory.
Finally, we conclude that this approach provides a nice tool to explain the different changes observed in the phase plane when changing a parameter (like in real systems) as a mixture of dynamic and geometric bifurcations.

Figure 2 .
Figure 2. One-plus-one geometric bifurcations exhibited in three-parameter spaces.Top: The bifurcation surface M reaches a maximum value at ε = 0 .Exploring the surface with horizontal planes, an equivalence class (cyan) disappears because the limiting bifurcation curve (red circle) collapses to a point when ε = 0 (top right panel).Bottom: The bifurcation surface has a saddle point.Passing through ε = 0 , the bifurcation curves (red branches of a hyperbola) change their position with respect to the asymptotes.As a consequence, the equivalence class around the central point switches (bottom right panel).

Figure 3 .
Figure 3. Two-plus-one geometric bifurcations exhibited in three-parameter spaces and subsidiary effects.In both cases C is a codimension-two bifurcation curve and M is a codimension-one bifurcation surface with a boundary along C. Left: Taking ε-slices two bifurcation curves (red) join in a unique one.Right: A bifurcation curve disappears after it collapses at the geometric bifurcation point.

Figure 4 .
Figure 4. (a) A one-plus-two geometric bifurcation.From the point (ε 1 , ε 2 ) = (0, 0) two one-plus-one geometric bifurcation curves (red) emerge.They split the parameter space into two regions corresponding to non-equivalent bifurcation diagrams.(b) Vertical one-dimensional "slices", obtained by fixing ε 1 and ε 2 , are used to explore the bifurcation diagram exhibited in the panel (a), that is, we consider the three-parameter space as a two-parameter family of bifurcation diagrams in a one-parameter space. https://doi.org/10.1038/s41598-024-61574-6

Figure 5 .
Figure 5. Theoretical scheme of the different topologies (Morse classification) of the bifurcation surfaces and the possible geometric and dynamic bifurcations on them.

Figure 6 .
Figure 6.Some extra possibilities of mixing different Morse surfaces and codimension-two bifurcation curves.

Figure 8
Figure 8. Schematic processes of some codimension-two-plus-one bifurcations in the HR model increasing the parameter ε .On the top, SC pictures with some bifurcations illustrate the global panorama at several values of ε .On bottom pictures some numerical bifurcation curves provide an outline of different ways, through "visible" vs. "invisible" folds, to arrive to similar global situations.Note that ε increases from left to right.

Figure 9 .
Figure 9. One-plus-one and one-plus-two homoclinic geometric bifurcations in the HR model.(a) (b, I) projections of the hom (1,2) bifurcation curve for different values of ε .The location of the codimension-one- plus-one bifurcation points ("visible folds") and the codimension-one-plus-two bifurcation point are given.(b) Schematic picture of the geometry of the hom (1,2) bifurcation surface presenting a cusp catastrophe geometry giving rise to the codimension-one-plus-two bifurcation point.(c) SC plots with the hom (1,2) bifurcation curves for three values of ε showing the change of their geometry.

Figure 10 .
Figure10.Three-parameter plot showing codimension-one homoclinic hom(1,2) bifurcation surfaces for the FitzHugh-Nagumo system.In plot (a) it is shown the theoretical structure of the homoclinic bifurcation surfaces, in this case it corresponds to the "duck-foot shape" case IIb of Fig.5.Plots (b) and (c) show the homoclinic bifurcation curves for several values of the small parameter ε and the codimension-two Belyakov bifurcation curve obtained using the continuation software AUTO.Plot (b) is given in the three-parameter space (α, s, ε) while plot (c) is represented by using the parameters (α, ε) and the AUTO L 2 -norm.